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Abstract. Dynamics of non-Brownian flexible fibers in Poiseuille flow between two parallel planar solid 
walls is evaluated from the Stokes equations, solved numerically by an accurate multipole code hydro- 
multipole. Fibers migrate towards a critical distance from the wall z c , which depends significantly on the 
fiber length N and bending stiffness A. Therefore, the calculated values of z c can be used to sort fibers. 
Three modes of the dynamics are found, depending on a shear-to-bending parameter F. In the first mode, 
stiff fibers deform only a little and accumulate close to the wall, as the result of a balance between the 
tendency to drift away from the channel and the repulsive hydrodynamic interaction with the wall. This 
mechanism is confirmed by simulations in the unbounded Poiseuille flow. In the second mode, flexible fibers 
deform significantly and accumulate far from the wall. In both modes, the tumbling pattern is repeatable. 
In the third mode, the fibers are even more curved, and their tumbling is irregular. 
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1 Introduction 

Dynamics of flexible fibers in simple shear and Poiseuille 
flows has been analyzed theoretically, numerically and ex- 
perimentally in numerous publications [Tl[2l[3ll4l[o1[6l[7l[8ll9l 
ITP1ITT1IT2] . Migration of fibers in Poiseuille flow [T^IT4TT5l 
16,l7j is the fundamental problem of modern lab-on-chip 
hydrodynamics, important in various biological, medical 
and industrial contexts, such as Brownian dynamics of 
proteins, actins, DNA or biological polymers, cell motion, 
swimming of microorganisms, drug delivery, transport of 
microparticles QMM15]. 

For significant pressure differences, corresponding to 
large maximal flow velocities, migration is caused by a 
fluid inertia |21) . However, fluid flows in microchannel de- 
vices often take place at low-Reynolds-numbers. In such 
systems, Brownian rigid rods migrate towards the wall [131 
[H] . and flexible fibers to an off-center position [22,23,24, 
[25]. 

For non-Brownian systems, the key question is un- 
der what conditions there exist off-center distances from 
microchannel walls where flexible fibers tend to accumu- 
late, what are their values, and how they depend on the 
fiber size, aspect ratio and flexibility. The importance of 
this problem is straightforward. Focusing of micro and 
nanoparticles is essential for their counting, detecting, and 
sorting [2"M2"7] . 

The dynamics of flexible fibers is also interesting from 
the fundamental point of view [35] . Evolution of their non- 



straight shapes is related to the existence of a family of 
modes, which are activated if the characteristic parameter 
exceeds subsequent threshold values. The parameter is de- 
termined as the ratio of the viscous forces to the bending 

ones [nEnunEa. 

In this paper, we study both practical and fundamen- 
tal aspects of the fiber dynamics. We investigate where 
the fibers accumulate, using the bead model and the mul- 
tipole method [33] of solving the Stokes equations, imple- 
mented in a very accurate, well-tested hydromultipole 
numerical code [30] . The goal is to determine how position 
of accumulation planes depends on the the fiber bending 
stiffness and its length, and to relate the findings to the 
characteristic parameter and its thresholds. In Sec. [2] we 
specify the system and theoretical model. The results are 
presented in Sec. [31 In Sec. 0] we conclude, discussing dif- 
ferent modes of the fiber dynamics and thresholds of the 
characteristic parameter. 



2 System 
2.1 Fluid flow 

We analyze motion and shape deformation of a single non- 
Brownian flexible fiber, moving freely in Poiseuille flow 
inside a channel made of two parallel solid walls, as illus- 
trated in Fig. [1] The fluid velocity v and pressure p satisfy 
the stationary Stokes equations [3T1I32] . 
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where ry is the fluid dynamic shear viscosity. 

The fluid is confined between two parallel infinite solid 
walls at z — and z — h, with the Poiseuille flow velocity 



F = F. e +F°, with 



v = Az(h - z)/h 2 St. 



(2) 



The stick boundary conditions are satisfied at the surface 
of the fiber and at the solid walls, which confine the fluid. 
At infinity, the fluid velocity v = vo. 

Distances are normalized by the fiber thickness d, ve- 
locities by the maximal velocity v m of the Poiseuille flow, 
forces by fo = 7H]dv m , and time by to = d/v m . 



wall 




centerline 



wall 



Fig. 1. A flexible fiber entrained by Poiseuille flow between 
two parallel solid walls. 



The system defined above is important for practical ap- 
plications, but complex to be studied theoretically. First, 
the shear rate depends on position z across the chan- 
nel, and second, the hydrodynamic interaction of the fiber 
with the walls is significant. To separate these two effects, 
we also study a reference system (see Fig. [2J, with the 
Poiseuille flow given by the same Eq. © , but not bounded 
by the walls, and extending beyond < z < h. 




centerline 



Fig. 2. Reference system: a flexible fiber entrained by un- 
bounded Poiseuille flow (without walls). 
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where k is the ratio of the Hooke's constant to fo and A 
is the ratio of the bending stiffness to fod 2 (in the follow- 
ing just called the bending stiffness). In the above equa- 
tion, Iq and U denote the equilibrium and time-dependent 
distances between the centers of the consecutive beads, 
respectively, with li = |tj|, where ti — ri—n-i is the dif- 
ference between the positions ru of the consecutive bead 
centers k = i— Here, ij = and V$ is the deriva- 
tive with respect to r,. The total non- hydrodynamic force 
applied to all the fiber beads vanishes, 5^ i=1 -Fi = 0. 

Translational and rotational velocities of the fiber beads, 
U = (Ui, Ujv) and = (Sli..., fijv), are linear combi- 
nations of the non-hydrodynamic forces F = (Fi, Fn) 
exerted on them all, and the multipoles of the ambient 
velocity field @, with the coefficients determined by the 
elements of the grand mobility matrix [34]. All the terms 
related to the ambient flow can be interpreted as result- 
ing from the hydrodynamic forces F = (Fox, F jy) and 
torques To = (Tqi, ...,Tbjv), exerted by the same ambient 
flow on motionless beads fixed at the same instanta- 
neous positions as the fiber beads, 



F + F 
T 



(5) 



with the mobility matrix (j, dependent on the instanta- 
neous positions of all the bead centers, r — (t\, rjy). 

For a given configuration, values of F , T and fj, are 
determined by the multipole expansion of the Stokes equa- 
tions [2{?ll3"5] , with the wall effects evaluated by the single- 
wall superposition [36 , 37] . The computations are performed 
with the use of the hydromultipole numerical code [30| . 
Then, the adaptive fourth-order Runge-Kutta method is 
applied to determine the fiber dynamics, 



dr/dt = U. 



(6) 



Initially, the fiber is aligned with the flow (i.e. along the 
x axis), with the bead centers located at = (i/o,0, zq), 
for i = 1, . . . , N . Owing to symmetry, the fiber moves in the 
xz plane. The computations are three-dimensional, and no 
deformation of the fiber out of the plane is observed. 



2.2 Fiber dynamics 



2.3 Parameters 



A single fiber consists of N solid spherical beads of di- 
ameter d equal to the fiber thickness [33]. Owing to non- 
hydrodynamic constraints, the beads do not move apart. 
There are no non-hydrodynamic torques, and the non- 
hydrodynamic force exerted on each bead i = 1, N by its 
neighbors is the sum of the elastic and bending forces [H] , 



In the numerical simulations, we have used single values of 
the bead diameter (length unit), the channel width h, the 
Hooke's constant k and the equilibrium distance between 
the consecutive beads Iq, 



h = 50, k = 80, l = 1.01. 



(7) 
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A large value of k and small gap size (Zo — 1) between the 
beads are chosen to model compact fibers which practi- 
cally do not change their length while bending. 

Three values of the fiber length N (in our units equal 
to the number of beads, or the aspect ratio) have been 
considered, with the corresponding fraction of the channel 
width, L = N/h, explicitly given in Table [TJ For clarity of 
presentation, we focus on discussing in details the results 
obtained for N = 10. 

Table 1. The fiber length L — N/h as a fraction of the channel 
width for the fiber aspect ratio N used in the simulations. 



N 


5 


10 


20 


L 


0.1 


0.2 


0.4 



Computations have been performed for a wide range 
of the initial fiber positions z$ across the channel. The val- 
ues of the bending stiffness A ranged from 0.01 < A < 2, 
and have been chosen to observe thresholds for different 
modes of the dynamics. It is known [^I fTUlfTTlfT^] that the 
transitions between C, S and W modes arc associated with 
specific values of a dimensionless parameter, equal to the 
ratio of the viscous forces (proportional to the local shear 
rate) to the bending ones. This parameter is widely used to 
characterize systems, which are far from interfaces. How- 
ever, it is clear that under confinement (as in the system 



considered in this work), there are additional wall effects 
which may influence thresholds of the fiber dynamics. In 
this paper, we are going to study these effects, by com- 
paring our system (Fig.Q} with the reference one (Fig. [2]). 
We use two basic parameters N and A to describe the 
fiber evolution. For N= 10, we evaluate a simple shear-to- 
bending dimensionless number mentioned above, 

r = {h/2 - z m )/A, (8) 

and analyze its critical values at the thresholds, and their 
dependence on the distance z m from the fiber center-of- 
mass to the closer wall. 



3 Results 

3.1 Lateral migration and accumulation planes 

In Ref. [38], dynamics of fibers in the same system was 
analyzed, focusing on the migration towards the central 
plane of the channel and its dependence on the fiber stiff- 
ness, aspect ratio and distance from the wall [38]. But for 
certain values of these parameters, fibers migrate away 
from the central plane. In this paper, we determine the 
critical distance z c from the wall where fibers accumulate. 

Fig. [3] shows evolution of the distance z m (t) between 
the fiber center-of-mass and the closer wall, starting from 
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Table 2. The distance z c from the wall where fibers accumulate. 



N\A 


0.025 


0.05 


0.125 


0.2 


0.25 


0.38 0.46 0.5 


1.0 


2.0 


5 


16.81 
3.6 


9.11 


4.48 


4.3 


4.19 


4.02 


4.3 


4.7 




(mi) 


(=bD.05) 


(=bD.02) 


(mi) 


(=t£).02) 


(=tD.05) 


(mi) 


(mi) 


10 




22.2 
10.1 


19.33 


14.45 


12.4 


8.9 8.3 8.2 


7.5 


6.4 






(mi) 


(mos) 


(m05) 


(J0.2) 


(mi) 


(mi) 


(mi) 


20 










21.28 

(mi) 


19.1 

(ifl.2) 


15.9 

(mi) 


14.2 

(mi) 



different values z m (0) = zq. Online, positions of the fibers, 
which move towards (away from) the central plane of the 
channel are plotted in red (blue). 

All fibers tend to an off-center position across the flow. 
For N = 5, the migration rate is significantly slower than 
for A^ = 20 (notice the f times larger range of times of 
the upper plots in Fig. [3]). 

The lateral migration of fibers is superimposed with 
oscillations of their center-of-mass position, clearly visible 
in Fig. [3] for A^ = 20. These oscillations are related to the 
fiber tumbling motion, caused by the local shear of the 
flow. A flipping time i/ is defined as the instant when the 
end-to-end vector (which links the centers of the first and 
the last beads of the fiber) becomes perpendicular to the 
flow direction. The distance from the fiber center-of-mass 
to the wall at time tf will be denoted as Zf, 

z f = z m {tf). (9) 

In Fig. [21 consecutive positions Zf are marked at four se- 
lected simulation runs with A" = 20. 

For all fibers, Zf tends to a critical position z c , 

Zf — > z c when t — > oo. (10) 

The value of z c depends on the fiber aspect ratio N, and 
the fiber stiffness A. Some of the evaluated values of z c 
are listed in Table [5J Their relative accuracy (typically, 
0.5-2%), is determined as the maximum of the oscillation 
amplitude and the half-a-distance between the closest de- 
creasing and increasing curves z m (tf), at the last flipping 
instant tf observed in our simulations. 

For A" = 10, the results are shown in Fig. [4] For 
larger values of A, the accumulation plane is located at 



N=10 

25 1 ■ ■ 




0.5 1 1.5 2 

A 

Fig. 4. Position z c of the accumulation plane versus the fiber 
bending stiffness A, for N — 10. 



the position z c larger than half of the fiber length, N/2, 
but smaller than A" (circles, violet online). However, for 
a smaller values of A, the accumulation distance rapidly 
increases with the decreasing A (squares, green online). 
When z c becomes sufficiently close to the mid-plane of 
the channel, the second accumulation plane is observed 
for the same value of A (diamond, magenta online). For 
N = 5 and A~ = 20, a similar tendency is visible in Tabled 

It seems that accumulation of stiff fibers is caused by 
the wall, which prevents them from escaping. Flexible 
fibers, however, accumulate far from the wall, probably 
owing to their shape deformation and the flow curvature. 
This hypothesis will be verified in the next section. 



Table 3. The position z"°" wa where fibers accumulate in the absence of walls. The arrows \ indicate that there is no 
accumulation position - all fibers migrate away towards z m < 0. 



N\A 


0.025 


0.05 


0.125 


0.2 


0.25 0.38 


0.46 
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2.0 
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(mi) 


8.65 

(mi) 


N 
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10 




22.2 
5.8 

(mi) 


17.85 

(mos) 


14.2 
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N=10 A=0.2 



N=10 A=2 




x 10 



N=10 A=0.2 





x 10 



N=10 A=2 



x 10 




x 10 



Fig. 5. Evolution of the distance z m (t) from a fiber center-of-mass to the closer plane, where the ambient Poiseuille flow 
vanishes. Top: the flow bounded by the walls, as in Fig. [T] Down: the same fiber and the flow, but without walls, as in Fig. [2] 



no - wa11 , see the left pan- 
However, the motion 



3.2 Comparison with unbounded Poiseuille flow helps 
to discriminate between two modes of accumulation 

In Table [21 we evaluate positions z^°~ wan of the accumula- 
tion planes for the same fibers and the same ambient flow 
as in the previous section, but without walls (the system 
shown in Fig. [5]). For more flexible fibers, the accumula- 
tion planes are located in approximately the same position 
with and without the walls, z c w z^ 
els of Fig. [5] (N = 10 and A = 0.2) 
of more stiff fibers significantly depends on the presence or 
absence of the walls. The difference can be seen by com- 
paring the top and bottom right panels of Fig. [5] (N = 10 
and A = 2). With walls, the fibers accumulate at Zf = z c 
inside the channel. Without walls, the fibers migrate out 
of the "channel regime" (defined as < Zf < h), what- 
ever is their initial position across the flow; there is no 
accumulation points in this range of Zf. In Table [21 such 
a behavior is indicated by arrows pointing down-right. 

Clearly, there exist two modes of the fiber accumula- 
tion inside the channel: caused by its hydrodynamic in- 
teraction with the wall (violet circles in Fig. [H F < 36) 
and caused by its interaction with the flow itself (green 
squares in Fig. [4l r > 42), where the parameter r is de- 
fined by Eq. ©, with the center-of-mass (always in this 
paper) taken at the flipping instant, z rn = zt. The tran- 



sition between both modes takes place for a critical value 
Jo such that 



36<r <42. 



(11) 



In the next section, we will investigate if the transition 
between both accumulation modes is correlated with a 
change of fiber shapes. 
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3.3 Fiber shapes 

In Fig. we compare evolution of fiber shapes for both 
modes. The snapshots are labeled by the corresponding 
values of the rescaled time, t/r, defined by the relations, 



t = t-t f {2), 
t = i/(3) - i/(2), 



(12) 
(13) 



where t/(n) is the instant of the n-th flip. The indicated 
fiber position Zf corresponds to the second flipping instant 
t f (2) (i.e. to i = 0). 

In the top panel of Fig. [Bl A = 0.2, and in the bottom 
one, A = 1. We first compare the snapshots taken at two 
critical positions z c from different accumulation modes. 
The first mode (accumulation caused be the walls, smaller 
r, violet color online), seen in the middle row of the top 
panel, corresponds to the S-shaped type of the motion. 
The second mode (accumulation independent of the walls, 
larger r, green color online), shown in the middle row of 
the bottom panel, is only slightly bended. 

To quantify this difference, we introduce two param- 
eters of a fiber shape: the curvature k (as in Ref. |11| ) 
and the fractional compression a (as in Ref. [12] )■ The 
time-dependent fiber curvature, 



N=10 




Fig. 7. The fiber curvature « at the instant of flipping at the 
accumulation plane z c , versus its bending stiffness A. 



N-l 



N 



r . ' 



2 n 

i=2 



(14) 



is defined as the mean inverse radius 1/r, of the circle 
determined by the centers of three consecutive beads. 

In Fig. [JJ we present values of k for fibers made of 
iV = 10 beads, with the center-of-mass at the accumula- 
tion plane z c at the flipping instant tf. 

The fiber fractional compression is defined as 12 



= 1 - 6(t f )/N, 



(15) 



where 6(tf) is the end-to-end distance of the fiber located 
at the accumulation distance z c in the time of flipping tf. 



(The end-to-end distance means the distance between the 
centers of the first and the last bead.) 

In Fig. [51 we present values of a for fibers made of 
N = 10 beads, with the center-of-mass at the accumula- 
tion plane z c at the flipping instant tf. 



N=10 




Fig. 8. The fiber fractional compression a, defined in Eq. (|15[l. 
at the flipping instant tf and distance z c , versus its bending 
stiffness A. 



From Figs. [JJ and [HI it is evident that for the second 
accumulation mode (A < 0.38, green squares), both the 
curvature and the fractional compression of the fibers flip- 
ping at the accumulation distance z c are much higher than 
for the first one (A > 0.46, violet circles), with the rapid 
change in the transition range, in agreement with the pre- 
vious analysis of the corresponding snapshots in Fig. [6] 

Until now, we have discussed the modes of the fiber 
dynamics only at the accumulation trajectories. In Figure 
[SI the fiber shapes at other trajectories are also shown. 
Eqs. (JHJ) and ([Til are now used to determine values of the 
parameter r for each trajectory. It is interesting that for 
r < r (the rows 1, 4, 5 and 6), the C-shaped type of 
motion is observed, and for r > Tq (the rows 2 and 3) - 
the S-shaped one. Such a transition to shape instability, 
triggered by a critical value of the shear-to-bending num- 
ber (equivalent to our To) is known in the literature, see 
9,12 and the references within. 



3.4 Third (irregular) mode of the fiber dynamics 

In the previous sections, only two modes of the fiber dy- 
namics have been discussed. However, in Figs. [JJ and [5J 
for a very small value of the bending stiffness A — 0.05, 
there appear also an accumulation plane (diamond, ma- 
genta online), which corresponds to much higher curva- 
ture and fractional conversion than the other ones. The 
corresponding (very compact) fiber shapes are shown in 
the lowest row of Fig. [9] From Fig. [4] it is clear that this 
plane is much closer to the wall than the accumulation 
planes of the second (green) type. Moreover, it is one of 
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Fig. 6. Evolution of a fiber shape (drawn to scale) for N = 10. Up: A = 0.2. Down: A — 1. Snapshots from simulations taken 
at the indicated times i/r (with i = t — tf(2) and r = t/(3) — t/(2), where i/(ra) is the instant of the n-th flip). The indicated 
values of Zf are attained at t = (second flip). 
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Fig. 9. Evolution of a fiber shape (drawn to scale) for AT = 10. Top: A = 0.2 and z c = 14.45. Middle: A = 0.05 and z ci = 22.2. 
Bottom: A = 0.05 and z C2 — 10.1. Snapshots from simulations are taken at the indicated values of the time i/r, defined by 
Eqs. (|12[1 - (|13| ). Zero and one correspond to the consecutive flipping instants. 
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two accumulation planes observed for the same value of 
the bending stiffness A = 0.05. The migration to these 
two planes is shown in Fig. 1101 with the accumulation dis- 
tances z c = 10.1 for the third, and z c — 22.2 for the second 
mode. 



N=10 A=0.05 
25r- ■ ■ ■ — 



E 

N 




2 4 6 8 10 

1 x10 4 

Fig. 10. Evolution of the distance z m (t) from the fiber center- 
of-mass to the wall, for N — 10 and A = 0.05. The lower 
accumulation plane corresponds to the third mode, and the 
upper one to the second mode. 



curve, with even larger value r « 900, is completely ir- 
regular. 



,■[--■ A=0.2 A=0.05 - - A=0.01| /V 




2 4 6 

i/r 

Fig. 11. The time-dependent curvature k of fibers with N = 10 
moving at z m . Dashed line (green online): z m fa z/ = 14.45 
and A = 0.2. Solid line (magenta online): z m ~ 2/ = 10.1 and 
A = 0.05. Dashed-dotted line (magenta online): z m ~ 16 and 
A = 0.01. The reduced time i/r is defined by Eqs. (fT2l - (fT3)l . 



The essential difference between trajectories correspond- 
ing to both modes is their time-dependence: regular for 
the second, and irregular for the third mode. This effect 
is visible in Fig. [TU] as small irregular fluctuations of the 
lower trajectories. This property can be used to determine 
the critical value -Ti of the parameter r at the transition 
between these two modes, 

260 < A < 298. (16) 

To study the nature of the irregular behavior, in Fig. [9] 
we compare evolution of shapes. At the trajectories of the 
second type (top and middle panels^ fibers bend in a re- 
peatable way, changing pattern almost periodically, with 
the half-period determined by the tumbling time r be- 
tween consecutive flips, which occur at integer values of 
i/r, c.f. Eqs. (|T2|) - (fT3|) for the notation. In contrast, the 
shapes of fibers at the third mode (bottom row) are not 
repeatable. 

In Fig. 111! we study time dependence of the fiber cur- 
vature. Colors visible online mark different modes of the 
fiber dynamics. For a fiber motion of the second mode 
(dashed curve, r « 53), k is almost perfectly periodic. The 
solid curve corresponds to the fiber shapes shown in the 
bottom row of Fig.[9l with r « 298 (the third mode), and 
it is quasi-periodic (but not regular). The dashed-dotted 

1 Notice that the fiber evolution shown in the top panel of 
Fig- EI with r w 53, is S-shaped, but in the middle panel, where 
F ~ 56, it is C-shaped. This example indicates that there is 
no universal correlation between the value of r and the shape 
type. This might be related to the value of z c = 22.2 very close 
to the middle plane of the channel, or other reasons. 



3.5 Tumbling time 

It is interesting to determine how tumbling of fibers de- 
pends on their position across the channel. In Fig. [TJ1 the 
fiber flipping frequency 1/r is plotted as a function of the 
distance Zf from the wall at the flipping instant. For the 
first and the second modes, it is a monotonically decreas- 
ing function (except fibers, which are close to the wall). 
Fibers at a larger distance from the wall tumble at a slower 
rate. Short fibers tumble more frequently in comparison 
to long fibers (notice a different scale on vertical axis of 
each panel in Fig. IT2"j) . Both effects are significant. The 
tumbling frequency is a bit larger for a smaller bending 
stiffness A. The difference is more pronounced for longer 
fibers. 

Following the idea of Brctherton [55], widely used in 
various contexts 1,2, we are now comparing the tum- 
bling time r characteristic for our flexible fibers entrained 
by the Poiseuille flow between two walls with the classic 
result of Jeffrey for the rotation half-period T/2 of rigid 
ellipsoids of revolution immersed in a simple shear flow in 
an unbounded fluid |40) . Jeffrey derived the following re- 
lation between the rotation frequency 2/T and the shear 
rate 7, 



T 7r(r e + l/r e )' 

where r e is the aspect ratio for the ellipsoid of revolution. 

For the Poiseuille flow, the shear rate depends on the 
position z across the channel, 

j(z) = 8(h/2-z)/h 2 , (18) 
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Fig. 12. The fiber tumbling frequency 1/r versus its distance 
Zf from the wall (solid lines), in comparison to the inverse half- 
period of a rigid spheroid with the aspect ratio TV (Jefferey, 
dashed line). Top: N = 5. Middle: TV = 10. Bottom: TV = 20. 



where, in our case, the channel width h = 50. 

In Fig. 1121 we compare our numerical results for the 
tumbling frequency of flexible fibers, 1/r, plotted as a 
function of Zf (solid lines), with the Jeffrey's linear rela- 
tion j(zf) /n(r e + l/r e ), which follows from Eqs. (fT7 |) -(fI5 |) . 
In general, an agreement would be expected for an effec- 



tive value of r e . As the reference, we plot the dashed line, 
which corresponds just to r e = TV. 

From Fig. [T2J it follows that for larger distances from 
the wall, 1/r is indeed proportional to (h/2 — Zf). The sur- 
prising effect is that the slope is quite well- approximated 
assuming that the effective hydrodynamic aspect ratio r e 
of flexible fibers is just equal to the number of beads TV, 

r E = N. (19) 

Unexpectedly, a better agreement is observed for longer 
and more flexible fibers which deform significantly during 
the tumbling, with the average geometrical aspect ratio 
much smaller than TV. 

For smaller distances from the wall, it is known from 
the literature that the Jeffrey approximation is not suffi- 
cient owing to the hydrodynamic interaction between the 
fiber and the wall, see Fig. 5 in Ref. [2J. 

The above discussion has been performed for the first 
and second mode of the dynamics (violet and green curves 
online). The third mode (magenta online) is seen in Fig. 1121 
as non-smooth, rapidly fluctuating lines, what reflects well 
the nature of this irregular mode. 



4 Conclusions 

In this work, we have considered dynamics of fibers, which 
are immersed in a low-Reynolds-number Poiseuille flow 
between two parallel planar solid walls at z — and z = h, 
and are initially aligned with the flow. Our key finding is 
that fibers with a different length (i.e. a different num- 
ber of segments TV) and a different ratio A of the bending 
stiffness to the flow amplitude, tend to accumulate at a 
different critical distance z c from the wall. The differences 
are pronounced. The dependence of z c on A and TV has 
been determined numerically in a wide range of the pa- 
rameters, based on more than 400 simulation runs. 

There exist two different mechanisms of the fiber accu- 
mulation. For stiff fibers, hydrodynamic interaction with 
the close wall prevents them from drifting out of the chan- 
nel. Therefore, in this case z c is a bit more than half of the 
fiber length TV, but still less than TV. This mechanism of 
the accumulation has been confirmed by the simulations 
performed for the same fiber and the same Poiseuille flow 
inside the channel range < z < h, but in the absence 
of walls. Without walls, stiff fibers migrate away from the 
channel range < z < h, whatever the initial position 
is. In contrast, flexible fibers tend to accumulate at larger- 
distances z c > TV, with similar values in the presence and 
in the absence of the walls. In this case, the accumulation 
mechanism is an interplay of the flow curvature, the fiber 
length and the fiber bending stiffness. 

The comparison with the unbounded flow was used as 
the criterion to discriminate between two different modes 
of the fiber dynamics. The additional differences between 
these two modes are the following. The first one has a 
larger fiber curvature and a larger fractional conversion. 
A correlation of the first mode with the C-shaped type 
of the motion, and of the second mode with the S-shaped 
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type of the motion has been observed for some cases. This 
problem will be studied in details elsewhere. 

A third mode has been also detected. Its basic feature 
(used to distinguish it from the second mode) is the ir- 
regular time dependence, best visible while analyzing the 
tumbling time, see Fig. [T^J In addition, for the third mode 
the fiber curvature and its fractional conversion are much 
larger than for the second mode. 

Our findings agree well with the previous literature re- 
lated to the fiber transitions to higher modes [ST ITUllTTlfT^l 
|2"B] . Following these papers, for N = 10 we have determined 
values of a shear-to-bending parameter _T = (h/2 — Zf) /A, 
and found its thresholds Jq and r% for the transitions be- 
tween the successive modes. Our findings are illustrated in 
Fig. 1131 Similar analysis has been (or can be) performed 
for other values of the fiber length N. This will be done 
in a next paper. 
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Fig. 13. Values of the shear-to-bending parameter _T for 
N=10. Three modes (different colors online) are separated by 
the thresholds Jo =36-42 and A =260-298. 



The results presented in this work can be used to 
sort non-Brownian flexible microfibers, depending on their 
length and bending stiffness. To this goal, additional mea- 
surement of their bending stiffness A is necessary. The 
analysis presented here indicates that neither the shape 
evolution nor the tumbling time is sufficient to determine 
specific value of A. 

Time and length scales of a fiber migration are rel- 
atively large. For example, in a microchannel of width 
h = 250 /im, with the maximal Poiseuille flow velocity 
v m —l mm /s, a fiber of thickness 5 ^m and length 100 /mi, 
initially located at the distance h/4 from the wall, typi- 
cally approaches a distance close to z c after 60 — 300 sec- 
onds, translating by 50 — 200 mm. 
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